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1 Introduction 



Inflationary cosmologyfl explains the oldness and the large-scale uniformity of the ob- 
servable universe. It may also explain the origin of the large-scale structure (galaxies, 
clusters of galaxies). It assumes the existence of a separate epoch in the cosmological 
history, when the expansion of the universe was extremely fast, "superluminal" . That 
epoch is called inflation. 

Inflation ends through a process known as reheating, in which the energy is trans- 
ferred from a homogeneous scalar field, the inflaton, to other fields and particles. In 
inflationary cosmology, all the matter in our universe has been created as a result of 
this process. 

Reheating can be viewed as amplification of quantum fluctuations of various fields 
in the time-dependent inflaton background. In slow-roll inflation, that background 
is an oscillating homogeneous inflaton field. It can be viewed alternatively as a Bose 
condensate of zero-momentum particles — the inflaton quanta. Then, reheating is the 
decay of the Bose condensate. 

It now appears probable that the postinflationary universe had been a much livelier 
place than it was thought to be. In many models of inflation, the decayof the inflaton 
condensate starts out as a rapid, explosive process, called preheating, a during which 
fluctuations of Bose fields coupled to the inflaton grow quasiexponentially and achieve 
large occupation numbers. This growth can be thought of as a result of parametric 
resonance. B S- H At the subsequent stage, called semiclassical thermalization,© the 
resonance smears out, and the fields reach a slowly evolving turbulent state with 
smooth power spectra. B B S @' 

The highly non-equilibrium states produced in the inflaton decay can support non- 
thermal symmetry restoration© and baryogenesis. Collisions of the classical 
waves produced in the decay process give rise to a strong, conceivably detectable 
background of relic gravitational waves. O 

Below, I discuss several topics in the theory of the resonant decay of the inflaton 
and in the description of the resulting non-equilibrium states. This discussion is 
based upon partially published^' S> work that I have done in collaboration with 
I. Tkachev. 



2 A model 

Consider the class of models in which the inflaton interacts with another scalar field 
X. The Lagrangian of the matter fields is 

L = \g» v d^d v <\> + \g^d,Xd v X - V{4>, X) (1) 



where 



X) = V+ffl + ^ 2 2 X 2 + l -M\X* (2) 
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and Vtj, is a potential for the field 0. The space-time metric is assumed to be that of 
a spatially flat Friedmann-Robertson- Walker (FRW) space-time with the scale factor 
a(t). In our simulations, the scale factor can be determined self-consistently, using 
the Einstein equations. 

For definiteness, let us now consider the so-called massive case, when 

^(0) = ^m 2 2 (3) 

Later, we will consider other potentials for 0. If we want our simple model to be 
a realistic model of inflation, the inflaton mass m should be of order 10 _6 Mpi, to 
comply with the COBE data. 

Both and X are regarded (for a while) as quantum operators, and we are 
working in the Heisenberg picture: the time dependence of averages is due to the 
time-dependence of operators, while the quantum state is time-independent. The 
system is translationally invariant in space. So, the average of the inflaton field in 
the quantum state is space-independent: 

(0(x,t)) = o (t) (4) 

We consider the initial data problem formulated as follows. 

The inflaton's zero-momentum (homogeneous) mode 0o rolls slowly during infla- 
tion, until it reaches a certain value, of order Mpi/3. At that time, the slow roll 
changes into oscillations of 0o about the minimum of V^. Thus the inflation ends and 
the reheating starts. The inflaton's homogeneous mode decays into fluctuations of 
X and those of its own field 0. The change from the slow-roll to the oscillatory be- 
havior of 0o is not instantaneous, and in principle the method described below allows 
us to follow the evolution of the fluctuations through the inflationary stage and into 
the reheating. The results, however, are practically unchanged if we adopt, instead, 
the following simpler formulation. 

To begin with, we solve for O and the scale factor a at the last stages of inflation 
and during the crossover from the slow roll to the oscillating stage. For this purpose, 
the fluctuations can be neglected. The crossover solution is insensitive to the values 
of 0o and d(f) /dt with which the slow roll starts. It turns out that the solution for 
0o, as a function of the conformal time rj, drj = dt/a(t), has an extremum during 
the crossover, before it makes the first oscillation. We will refer to the extremum 
as the time when the inflaton oscillations start; we take that time as rj = 0. In 
the massive model, 0o(O) = 0.28Mpi.@ The fluctuations at rj = are small, and in 
our formulation we assume that they are in the conformal vacuum, with respect to 
annihilation operators defined below. This specifies the (quantum) initial conditions 
for the fluctuations. (The assumption of the conformal vacuum is not a crucial one, 
as long as it gives initial fluctuations of the correct order of magnitude. @) 

It is convenient to rescale to the dimensionless conformal variables: 
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T = TflTj j (y9 = 0a/0o(O) 

£ = mx; x = Xa/0 o (O) 

We have normalized the scale factor a(r) to 1 at r = 0, when the inflaton oscillations 
start: a(0) = 1. By virtue of our rescaling, <£>o(0) = 1- The action for the matter 
fields in the new variables takes the form 



s $(0) 



J d^dr 



1 1 A/f 2 

1 2 2 o 2 2 1 iw X 2 2 

--a tp -2qip X ~ X 

2 2 m z 



(5) 



Hereinafter dots denote derivatives with respect to r; h = a/a is the rescaled Hubble 
parameter. 

The parameter q appearing in fl5|) is defined as 

and is called the resonance parameter. 

In this type of inflationary models, the resonance parameter q is naturally large. @ 
Indeed, with 0q(O) ~ O.lMpj and m 2 ~ 10 — 12 ikfp l5 one gets q ~ 10 10 g 2 , so even for 
fairly small g 2 one still has g> 1. A large value of q is necessary for efficient resonant 
decay of the inflaton in this model. @ 

The Heisenberg equations of motion following from (||) are 

(p - Vlp --tp + a 2 tp + Aqtpx 2 = (7) 
s a 

X ~ V|x - ^x + m 2 x a 2 x + M<p 2 X = (8) 
where m x = Mx / m. We can write 

^r)=^(r) + ^,r) (9) 

where 6ip(£, r) is the (quantum) fluctuation. Recall that 4>q{t) defined in (f|) and, 
consequently, <po(r) are c-number functions. 



3 Quantities to compute 

What are the interesting quantities to compute? One is the power spectrum of fluc- 
tuations. This is defined as follows. Introduce the spatial Fourier components 

X ^ T) = I jB^ Xk(T)e ^ ' Me T)= I (2% ^ {T)e ^ (10) 
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Thus k is a comoving momentum (in the units of the mass m); it does not redshift in 
an expanding universe. The power spectra are then defined as 



Px(k) = (\Xk\ )/V > P M = (\PkT)/V (11) 

where V is a normalization spatial volume. Because the system is isotropic, the power 
spectra are functions only of k = |k|. 

The angular brackets so far denote averaging over the quantum state of the system. 
The actual calculations were done using the classical approximation; we will come to 
this important point shortly. In the classical approximation, the averages over the 
quantum state are replaced by averages over realizations of the random initial data.0 

Another interesting quantity is the variance of a field. With the above definitions, 
it would be just an integral over d 3 k of the field's power spectrum, if we did not have 
to subtract away the high-momentum quantum fluctuations in order to make the 
variance finite. A natural scheme is to subtract the full initial value of the variance: 

( X 2 )i(r)= J[P x (k,r)~P x (k,0)}d 3 k (12) 

If we introduce a maximal (cutoff) momentum fc max , we can view ( |T2"D as the limit of 

( X 2 ) 2 (r)=/ P x (k,r)d 3 k- [ P x (k,0)d 3 k (13) 

as /c max goes to infinity. A maximal momentum is automatically present in our lattice 
calculations. We have to make sure that it is large enough so that ( jl3|) is insensitive to 
the cutoff or, equivalently, the modes with k ~ fc max do not give a large contribution 
to (0). 

The variances measure the typical size of the fluctuations. For example, they 
determine the fluctuations' symmetry restoring power. 110 Of main interest are the 
values of the variances at sufficiently large times, when the fluctuations are the largest. 
For the values of k max that we use, the first term on the right-hand side of (0), when 
resonance is efficient, quickly becomes much larger than the second. So, we will 
"forget" the second term and use 



( X 2 )(r)= P x {k,r)d 3 k (14) 

The parameters of our model (m, Mx, g 2 ) can be viewed as the running parameters 
of quantum field theory, normalized at the momentum point /c max . 

In our actual calculations, the variances are obtained as lattice averages: 

1 N s 

(X 2 ) = W Y,X 2 (15) 

ly s i=l 

«M 2 > = (16) 
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where i labels the lattice sites, and N s is their total number; 

I n s 

jj-J^<Pi = <Po (17) 

s i=i 

is the inflaton's zero-momentum mode. 

Note that in the classical approximation, the variances (|i~5|), (|i~6|) and the zero 
mode (|l7l) need not be averaged over realizations of the initial data. These quantities 
are self-averaging, as they are defined as averages over a large number lattice sites. 
Similarly, for the power spectra in the classical approximation, we use, instead of 

1 r , , dQ, , , 1 r , , dQ, , . 

AM = F /l*klV p ^) = f/W^ (18) 

where f2 is the direction of k. 



4 The linear stage 

As seen from (|5|), the parameter m 2 /0Q(O) ~ 10 -11 suppresses the initial magnitude 
of the fluctuations relative to the scale at which their non-linear interactions become 
important. Unless the parameter q is too large, in which case the scale of non-linearity 
is itself too small, the initial evolution of the fluctuations proceeds in the linear regime. 

Expanding the equations (0)~© to the first order in 5(p and \i we obtain the 
following linearized equations 

<fo- (a/a)ip Q + aVo = (19) 
5(p - (d/a)5cp + a 2 5ip = (20) 
X - V?X - (a/a) X + m\a\ + 4q<p 2 oX = (21) 

m x = Mx/rn. The Fourier expansion (0) turns the equation fl2lf ) into 

X k + ^ 2 Wx k = (22) 

where 

luI(t) = m 2 x a 2 (r) + k 2 - d/a + 4q<p 2 (r) (23) 

Eq. ( ^2|) is still a quantum equation, but it is a linear equation and can be solved 
exactly. 

The operator solution to (|22|) can be written as 

Xk(r) = fk(r)b k (0) + /;(r)6^ k (0) (24) 
The c-number function /^(t) satisfies the same equation as x\^i 

f k + u 2 k fk = (25) 
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and the initial conditions 

fk(0) = i — (26) 

0(0)^(0) 

/ fc (0) = [-iu k (0)+h{0)]f k {0) (27) 

where 

ul(r) = m\a\r) + k 2 + Aq V 2 (r) = u 2 k (r) + al/a (28) 

Note that fa depends only on the absolute value of k. The operators b and are 
time-independent creation and annihilation operators normalized according to 

[b k , b^} = <5(k - k') (29) 

The quantum state of the system is assumed to be the vacuum of b-^. The time- 
dependence of the operators X]c * s entirely due to that of the c-number mode functions 

fa- 

5 Non-adiabatic amplification 

In the linear approximation, the evolution of the inflaton zero-momentum mode is 
given by ([II]). In a static universe, i.e. for a(r) = 1, the relevant solution to (|19|) is 



<Po{t) = cos r (30) 

Then, the equation (p5|) reduces to the well-known Mathieu equation: 

/ fc + L4 + 2gcos2T]/ fc = (31) 

where A = k 2 + vr? + 2q. The instability zones of the Mathieu equation, in the plane 
of A and q, determine momenta k for which there are exponentially growing solutions. 
This growth is known as parametric resonance. The case of large q is known in the 
literature as the broad resonance case.B 

In the expanding universe, solutions to flT9"| ) are not periodic. The crucial point 
is that, when q is sufficiently large, there still are growing solutions (although the 
growth is no longer simply exponential), see Fig. [I]. The range of q for which there 
is an efficient growth of the variance of X has been determined numerically. @ For 
example, for Mx = 0, it is q 10 3 . The growing solutions in the expanding universe 
have been recently explored analytically. 

To distinguish the case with a non-periodic (fo from the case of a periodic (fo, we 
will use, instead of the term parametric resonance, the term non-adiabatic amplifica- 
tion. This stresses that the periodicity of <po is not crucial for the existence of the 
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Figure 1: Variances of the physical fields X and 0, together with the inflaton's zero- 
momentum mode, in the massive model with q = 10 4 , Mx = in the expanding 
universe. 



growing solutions. What is crucial for the efficient amplification is that, at least once 
in a while, the adiabatic (WKB) condition for Eq. (|22|) 

^ (T)I «1 (32) 



w2(r) 

breaks down. 

Soon after the inflaton oscillations start, the inflaton zero-momentum mode is well 
approximated by 

¥>o(t) ~ { p( T ) cos(ar + a) (33) 

where (p is a relatively slowly decreasing amplitude. During the linear stage, its 
decrease is entirely due to the redshift of the field: 

Cp{r) ~ a- 1/2 (r) (34) 

The amplitude of the physical field 0o is 4>{j) = 0o(O)<^(r)/a(r). 
It is convenient to introduce the effective resonance parameter!! 

q ^ T) = = q ^T) (35) 

For example, in Fig. Ill, at r ~ 10, we have 2 ~ 10 _5 Mp[, so q e g ~ 1. 

The most interesting case is that when the resonance parameter q is large enough 
for the amplification to fully develop. That means that the growth of the amplified 
modes is cut-off by non-linear effects, rather than by the expansion of the universe. 
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This is the case when the fluctuations can reach the largest values. It occurs when 
at the end of the linear stage the effective resonance parameter is still fairly large: 
q c g ?£ max{l, (Mjf/m) 4 }, which translates into q ?Z 10 4 max{l, (M x /m) 4 }M 

As long as q^ ^> 1, the adiabatic approximation is good most of the time. The 
condition (|3"2"|) breaks down only during small intervals when ip is close to zero: 



cos(ar + a) < g cfT 1/4 (36) 

and, even then, only provided that 

k 2 /a 2 + m 2 . < q^ 2 (37) 

These intervals of time correspond to spikes in (X 2 ), which are seen in Fig. |l|. We 
will distinguish between (x 2 ) in spikes and in "valleys" between them by a subscript: 
"s" or "v". One can show analytically that during the linear stage (x 2 )s/ (x 2 )v ~ ql^- 
The band of momenta specified by (1371) is analogous to the strongest (for a given 



q) instability band of the Mathieu equation (31). Notice that k/a is the physical 



momentum in the units of m. For a given k, even one that satisfies (p7J), the mode 
function will not grow each time the adiabatic condition is brokenJH' so the 
power spectrum of x a t the linear stage becomes jagged. 



6 Quantum-to-classical transition 



DC- 



The modes of x with growing achieve large occupation numbers and begin to 
have classically. The theory of this phenomenon is parallel to the linear theoryCJ' 
describing how large-scale classical density perturbations are generated from the vac- 
uum fluctuations during inflation. 

The occupation number in the modes with |k| = k, obtained by a Bogoliubov 
transformation, is 



(0) 



n k (r) = 2 jy {r) [\fk ~ hf k \\r)+ui(r)\f k \ 2 (r)\ - 1/2 (38) 

One can show that up to 0(l/rik) the quantum averages of the operators in the modes 
with |k| = k can be approximated by classical averages computed with the help of a 
certain distribution function. When the state vector is the vacuum of 6^ and 6_j £ , 

the distribution function for the pair of modes with momenta k and — k isJUi 

F\Xk, X k ; r] = M exp (- y^l^i 2 ) 6 (AWxk - Mr)x k ) (39) 

where the delta function is a shorthand for the product of the delta functions for real 
and imaginary parts; M is a time-independent normalization. (Note that half of the 
values of k give all the independent X\c-) 
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It is useful to notice that (|3TJ) satisfies the classical Liouville equation during the 
entire linear stage. So, we can evolve it back and use it as an initial condition at 
r = 0, where /& and its derivative are given by (p6[), (j2T|). This leads to the classical 
problem with random initial conditions distributed according to@ 

•Ffckk' *k; °] = X' ex P (-^^^(°)lxkl 2 J 6 (*k + " ^(°)]Xk) ( 4 °) 

V = K/(27r) 3 . This classical problem and the original quantum problem "converge" 
as the modes are amplified. 

How good is this classical approximation? For example, if q e g ~ 1 at the end of 
the linear stage, the typical occupation number in the efficiently amplified modes at 
that time is n amp ~ cf) 2 jm 2 . The corrections to the classical approximation are then 
of order l/n amp ~ 10~ 7 . 

When 

\fk\\fk\ ~m 2 n fc /0o(O) (41) 

and n k ^> 1, one can derive the classical approximation in the following simple way. 
The solutions /fc(r) and flij) to (p5|)-(|27|) satisfy 

h(r)r k (r) - r k (r)f k (r) = (m 2 /^(0))^ (42) 

Given (II]) and n k 3> 1, the right-hand side of ( f4"2"D can be neglected compared to 
each of the terms on the left-hand side. That means that the phase 9 k defined by 

f* k =f k eMWk) (43) 

becomes almost time-independent. With the definition (fl3"D, the operator solution 
( |24"D takes the form 

Xkto = [ b k + eldkb -k\ = (44) 

When the phase is almost time- independent, the time-derivative of %k can ^ e 
found as 

*k' ~ / fe '°k' (4 5 ) 
provided is large enough, cf. ( p]) . The canonical momentum conjugate to Xk> as 
obtained from the action (|), is proportional to (x_k — ^X_k)- A direct calculation 
shows that [(D(k), C?(k')] = 0, for any k and k'. So, in this approximation, Xk an d 
its canonical momentum commute and can therefore be regarded as random classical 
variables. 

Even if at some time the condition (O) breaks down, i.e. one of the two quantities, 
|/fc| or \fk\i happens to be small, the classical approximation will still work, as long 
as «fc > 1. The reason is that the other of the two quantities is necessarily large, cf. 
(p8|), so the quantum "noise" will not be able to disturb the evolution of the modes. 
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In the vacuum of 5^., the variance (|12|) is 




As long as q e g ^> 1, fj. oscillates rapidly in the "valleys". When Q k becomes almost 
time-independent, |/fc| 2 also oscillates rapidly there. As a result, the variance acquires 
a high-frequency "beard" , see Fig. [TJ. 



7 Rescattering 

The modes outside the band (R7|) and fluctuations of the field if do not undergo an 
efficient non-adiabatic amplification in the massive model. The system thus arrives 
into the regime when the low-momentum fluctuations of x are already large, while 
its higher-momentum fluctuations and the fluctuations of ip are still small. 

Fluctuations of <p now get amplified via rescattering® of the already large fluc- 
tuations in the low-momentum modes of X- Further rescattering of the x an d V 9 
fluctuations leads to amplification of the higher-momentum modes of both ip and 

x .0 

The early stages of these processes can be studied using the following approxima- 
tion. We assume that the fluctuations of ip are much smaller than the amplitude of 
its zero-momentum mode: 

<(M 2 > « f (47) 

but we make no assumptions about the size of the fluctuations of x- Expanding the 
full equations (0), (§|) in Sip (but no longer in x), we obtain, to the first order, the 
following approximate equations (p ^ 0) 



= (4* 



Xk + ^fc( r )Xk + 8 <?^o(r) j d 3 p<ppXk+p ~ 

<Pp + n 2 p (r)cp p + Aq<p (r) Jd 3 kxlxy: + p = (49) 
ip> -(d/a)(p + a 2 (T)ip + Aq(x 2 )(p = 



(50) 



where Q 2 = p 2 — a/a + a 

The main contribution to the last term on the left-hand side of ( 49|) comes from the 
low-momentum, classical modes of x- The fluctuations of ip are driven by a classical 
force and therefore become classical themselves. This quantum-to-classical transition 
happens in the linear regime with respect to 5<p. In this respect, it is similar to the one 
for the low-momentum modes of Xi even though it is due to the force-driven rather 
than parametric amplification. We obtain the initial distribution function for ipp and 
0p analogous to (0). In fact, within the accuracy of the classical approximation, we 
could have taken the initial values for ip and <p to be zeroes: the classical force would 
produce fluctuations of tp anyway. 
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The amplified fluctuations of ip, together with the low-momentum fluctuations of 
X form a classical force for the higher- momentum fluctuations of Xi see (EH)- m this 
way, the power spectra propagate to larger momenta. 



According to fl49|) , the 5<p subsystem is a resonator of high quality (the universe 
already expands slowly at this stage). Therefore, the rapid growth of 5tp will continue 
until 5(p ~ (p, and the approximation (f£8|)-(f£9D, based on (^7|), breaks down. A 
more precise condition for the end of the rapid growth of 5ip, obtained from the 
simulations, is 

<(M 2 > ~ o.i^ 2 (51) 

The system then enters the fully non-linear stage. 



8 The fully non-linear calculation 

To include all non-linear effects, we simulate the complete non-linear equations of 
motion (|7]), (|8|) as a classical non-linear problem with random initial conditions for 
X and Sip. The random initial conditions are distributed according to (f40|) and a 
similar expression for 8tp. The initial conditions for (p are: (po(0) = 1, 0o(O) = 0. 
The simulations were done on 128 3 lattices, with m 2 = 10 Mp 1; both for massless 
and for massive For massless X, we determined the evolution of the scale factor 
a(r) self-consistently, taking into account the effect of the produced fluctuations in 
the Einstein equations. For massive X, we used the matter-dominated form of a(r). 
Simulations were done also for the conformally-invariant model, discussed below. 
Additional results of the simulations are given in Figs. pi 0, H. The variances and 



the power spectra are defined according to (15), 

We used the classical distribution (|I]) for all momenta present on the lattice. This 
is possible because, even at the largest times that we consider, the main contributors 
to the physical quantities are modes with large occupation numbers. In general, 
a lattice calculation of the power spectrum can be trusted not up to the maximal 
momentum fc max present on the lattice, but only up to the Nyquist momentum k^ y ; on 
a cubic lattice k^y = & max / Vs. For example, in Fig. |2| (obtained for the conformally- 
invariant model), fc max ~ 28 and &N y ~ 16. For 3 < k < 16, the latest outputs of the 
power spectrum are close to a perfect straight line in the log-linear plot: 

P x (k) ^ P ex P (-k/k ) (52) 

Such exponential fronts are obtained also for the power spectrum of <p and the fields 
in other models. @' @ Integrals of modest powers of k computed with the power 
spectrum ( [5*2]) will be saturated at k ~ k , where P x {k) ~ P , and the occupation 
numbers are large. 

This distinguishes the highly non-equilibrium states emerging after the inflaton's 
resonant decay from the equilibrium, thermal states. For a real scalar field with a 
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Figure 2: Power spectrum of the field x hi t ne conformally invariant model, output 
every half-period at the extrema of <^o( r )- 



time-independent dispersion law e(k), the occupation numbers are given by the Bose 
distribution 7Zb(/c), and the power spectrum, at temperature T, is 

P T (k) = (2n)- 3 e-\k)n B (k) = ^'Y'^ ^ 

1 — exp[e(k)/ 1 \ 

The variance subtracted at T = is 

(x 2 )t - (x\ = f d'kP T {k) (54) 



The modes with e(k) <C T are essentially classical but the main contributors to fl54 ) 
are the quantal modes, those with tib ~ 1- 

The new states are to a good accuracy classical, precisely because they are highly 
non-thermal, with the power spectra behaving like (0). For these states, there is no 
ultraviolet catastrophe (no Rayleigh- Jeans problem). Eventually, the system will ther- 
malize, and will be dominated by quantum modes. Nevertheless, there is a prolonged 
stage in the evolution when the power spectra are concentrated at small momenta, 
and the system is essentially classical. 

The non-linear phenomena described by the classical equations of motion are cre- 
ation, scattering, decay, and annihilation of the classical waves arising through the 
amplification of quantum fluctuations. We refer to these phenomena as rescattering.lt 
In the language of particle physics, rescattering can be viewed as stimulated creation, 
scattering, decay and annihilation of "particles" in states with large occupation num- 
bers. A good analogy would be the scattering of waves from throwing pebbles in 
a pond, if those waves could be made high enough for their non-linear interactions 
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(rather than merely interference) to become important. Such non-linear interactions 
of waves occur in plasma physics and fluid dynamics. El' 

9 Approximate methods for the non-linear stage 

9.1 The Hartree approximation 

We will see below that as long as 

<(M 2 > « f/^ (55) 

we can neglect the fluctuations of if altogether. For large q e s, (|55|) is a stronger 
condition than ([IT]). In terms of the physical field <ft, ( |5o] ) reads 

<W) 2 >«0o(O)/<7 (56) 
We then get a system of two equations 

Xk + ^( r )*k = ( 5? ) 

y3 - (d/a)ip + a 2 (r)if + Aq( X 2 )^ = (58) 

These comprise the Hartree approximation for the present model. @> B The classical 
average in (^) approximates the corresponding quantum averageS with the accuracy 
0(l/n amp ), where n amp is the typical occupation number in the amplified modes. 



The equation (|57| ) is formally identical to fl22]) , so the operator solution (|24]) still 



applies, although the form of cpo and, consequently, of fk will now be different. The 



phase 9k, see (|3|), becomes almost time- independent under the same conditions ([41]) 
and rik ^> 1. So, as long as q e g 3> 1, the variance of \ in the Hartree approximation 
has the high-frequency "beard". 

The time r sc when the "beard" disappears on our pictures while q c s ^> 1 (e.g. 
t pa 10.5 in Fig. [3]) is the time when the Hartree approximation breaks down, and 
rescattering becomes a dominant effect. 

9.2 Corrections to the Hartree approximation 

To derive the condition (|55|), we need to estimate corrections to the Hartree approxi- 
mation and see when they become large. While ([|7]) still applies, we can use for that 
purpose the approximate equations (|4l|)-(|4"9"|). We will do the estimate for the general 
case of a field x with N real components. The components will be labeled by Greek 
indices a, (3 r . . ; a repeated index is summed over and is sometimes omitted. 
Eq. © is solved by (p ^ 0) 

Mr) = V%\t) + 4g £ dr'G p (r, t> (t') / d 3 k X * ka X k+p Jr') (59) 
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Figure 3: Same as in Fig. 1 but for q = 10 6 , Mx = 0. 



where ipp is the solution to the free equation (|20|) , and G p (t,t') is the retarded 
Green function for (|20|): G p (t,t') — fi^ x (r) sin[cUp(r — r')], where ^ ps + a 2 . 

The magnitude of <^p is of the order of the initial condition for (p-p and is very small. 
At the stage when the fluctuations of x still rapidly grow, the main contribution to 
the time integral in (|59|) comes from an interval of duration of order a~ l near r' = r, 
when the x fluctuations are the largest. We then estimate (|59|) as 

qCp{r) 



T 



(60) 



where we have used the estimated p ~ a(r) for the typical momentum of the ip 
fluctuations, so that also fl p ~ a(r). 
Eq. (PD is solved by 



Xkfs( T ) = X^( T ) + 8q J dT'F k {TiT')ip {T') J d 3 pifpXk + p^ 



(61) 



where xj^( r ) solves Eq. (|57|) of the Hartree problem 
retarded Green function for (1571); Fdr = r') = 0, Fi.fr = r' 



-(HD, and F k (r,T') is the 
—1. The integral over 



time in (|6T| ) is estimated in the same manner as the one in (|59| ) 



awfc(r) 



d WpXk + p„ 



(62) 



where u k (r) ~ y/q<f{r). 

Using the estimates (|6~0"D and (|62l) in ( |6lf ) we obtain 



(o) gV 



rf 3 fc^ 3 PXk' + p,^k',^k + p^ 



(63) 
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We recall at this point that x an d <fi are random variables. By the translational 
invariance and the O(N) symmetry, 



(XkA'p) = ^(XkXk^kk' "> (2vr) 3 iV-V- 1 (x^k)^^(k - k') (64) 

So, for instance, the average of tp-p (p 7^ 0) over realizations of the initial data is zero. 

The first correction to the Hartree approximation is obtained when we replace 
X on the right-hand sides of floPf), ( 161]) with x ■ m the Hartree approximation, 
different pairs of modes of x with the opposite momenta (k and — k) are statistically 
independent. For example, for any p / 

/ *(o) *(o) (0) (0) v = 

\ A k/3 A k'+p,a A k'a' A k+p, / 3' / 

/ *(o) *(o) w (o) (o) v / *(o) (0) w *(o) (o) v r6 ,x 

for any p / 0. Higher-order correlators decompose in a similar way. 

Using this approximation in (j60|) , we estimate the mean square of \ipp\ as 

(i ^p |2) ~ Swl d3k{ix ^ 2){ix ^p i2) (66) 

and the variance ((8(p) 2 ) = V~ x J d 3 p\{pp\ 2 as 

((M 2 )~?V(x 2 )> 4 iV (67) 

We can similarly estimate the first correction to (xjLXk)- Towards the end of the 
Hartree stage, it is of order of the average square of the right-hand side of (|63|) . (The 
cross correlation between the right-hand side and at these times a smaller effect 

for large iVj We find that the first correction becomes of the order of the leading 
term when lit 

According to (|6"TD, at that time 



((M 2 ) ~ ^ i-e. (0 2 )~^(O)/g (69) 

cf. the condition (|56|). 

The estimate (|68| ) grows indefinitely with AT, for a fixed g. The physical variance 
cannot do that. So, for sufficiently large N, N > N Q (q), where N Q increases with 
q, the Hartree approximation does not break down at all during the rapid growth of 
(x 2 ), i-e. the condition (|69|) is never reached. 

Otherwise, the estimate fl6"P| ) determines the time r sc after which rescattering be- 
gins to strongly influence the evolution of (x 2 )- Because in deriving ( |S9|) we neglected 
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various numerical factors, it should be more of a scaling law with respect to q than 
of an actual numerical estimate. Nevertheless, it turns out to correspond numeri- 
cally quite well to ((5ip) 2 ) at the time when the "beard" disappears on our pictures, 
whenever that happens while q e & ^> 1. For example, in Fig. ||, the beard disappears 
at r = r sc « 10.5. At that time, ((5<fi) 2 ) ~ lO^M^, in agreement with On 
the other hand, in Fig. [I], q e s becomes of order 1 already at r ~ 10, and even the 
maximal value of ((5<f)) 2 ) is somewhat smaller than the estimate flSDp. 

At t > r sc , the system exhibits chaotic behavior, characteristic of a non-linear 
classical system. The stage when 

^(O)/g<((50) 2 )«O.10 2 (70) 

has been called the chaotization stage. 0' S When the second condition in ([70]) is 
broken, the system enters the fully non-linear stage with developed chaotic behavior. 
After chaotization, the system is in a slowly evolving (quasi-steady) state with smooth 
power spectra, see e.g. Fig. 0. Such states are naturally called turbulent. It is the 
chaotic evolution that finally establishes the maximal values attained by the variances 
of the fields. 

The turbulent states emerge in a variety of models of the inflaton's resonant 
decay (cf. the conformally invariant case below). The power spectra in these states 
in general consist of a segment of a power law at the smallest k, k ~ 1, from where 
the exponential front ( |52"D takes off. We expect that at later times, the front will 
slowly move to larger momenta, and the region of k occupied by the power law will 
grow. 



10 Decay of the inflaton 

An interesting feature of the massive model in the expanding universe is that the 
inflaton decays completely fairly soon after it starts to oscillate, see Figs. II |3|. 
This early decay does not occur in the same model without the expansion or in the 
conformally invariant case (discussed below). 

Because the inflaton decays fast, and the expansion of the universe is already 
slow, the total energy is almost conserved during the decay. One can obtain a simple, 
albeit non-rigorous, estimate of the maximal size of X fluctuations from the energy 
conservation. We have (in terms of the physical fields) 

m^fa) ~ gPiiStWfaXX 2 )^ (71) 

where T\ is a moment right before the decay, and r 2 is the moment when (X 2 ) is 
maximal; we have assumed that a significant portion of the inflaton's energy went 
into the X fluctuations. Using ((<50) 2 (ti)) ~ 0.1(/> 2 (ti), see (|51|), and assuming that 
((5(p) 2 ) does not change drastically during the decay, we get 

<(ty) 2 >(7a) "0.1^(71) (72) 
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From ([Tl|) and (|72|), (X 2 ) max ~ m 2 / g 2 ~ 0Q(O)/g (neglecting numerical factors), or 

(X 2 ) max = e<Pl(0)/q (73) 

where e takes into account the fraction of the inflaton's energy that went to X and 
other numerical factors. 

The 1/q estimate ( |T3"D agrees with the scaling of our data points and with an 
estimate obtained from the self-similarity hypothesis. From the data, we find e ~ 
0.1-1. 

11 Conformally- invariant model 

This model has the scalar potential 

y(0,X) = ^ 4 + ^</> 2 X 2 (74) 

The fields 4> and X are massless. For the inflaton self-coupling A we take a realistic 
value A = 10" 13 . 

The inflaton oscillations start at 0o(O) ~ 0.35Mpi , d(po(0)/dr] = 0;^ r] is the 
conformal time: dr] = dt/a(t). The rescaled conformal variables are 

r = a/A0 o (O)?7 ; tp = 0a/0 o (O) 

£ = v / A0o(O)x ; X = Xa/MQ) 

so that y?o(0) = 1. The universe is radiation dominated, a(r)/a(0) ~ 0.51r + 1. 
The equations of motion in the conformal variables 

(p - V\y + if* + 4g X V = , 

X - V|x + Wx = • (75) 

are the same as in the flat space-time. The resonance parameter in this case is 
q = g 2 /AX. 

The power spectrum of x f° r Q — 30 is shown in Fig. |[ Notice the formation of 
the exponential tail, cf. (|52|) . The latest outputs in Fig. |2| fall almost exactly onto 
each other — a sign of a quasi-steady state. 

12 Conclusions 

Resonant decay of the inflaton leads to formation of highly non-equilibrium states 
after inflation. These states are dominated by interacting classical waves. As we 
have shown elsewhere,® collisions of these classical waves give rise to a potentially 
observable background of relic gravitational waves. 



18 



Acknowledgments 



I would like to thank the organizers of SEWM'97 for holding a stimulating conference. 
The author's work was supported in part by the U.S. Department of Energy under 
Grant DE-FG02-91ER40681 (Task B), by the National Science Foundation under 
Grant PHY 95-01458, and by the Alfred P. Sloan Foundation. 

References 

[1] A. D. Linde, Particle Physics and Inflationary Cosmology (Harwood, Chur, 
Switzerland, 1990); E. W. Kolb and M. S. Turner, The Early Universe (Addison- 
Wesley, Redwood City, California, 1990). 

[2] L. Kofman, A. Linde, and A. A. Starobinsky, Phys. Rev. Lett. 73, 3195 (1994). 

[3] J. H. Traschen and R. H. Brandenberger, Phys. Rev. D 42, 2491 (1990); 
Y. Shtanov, J. Traschen, and R. Brandenberger, Phys. Rev. D 51, 5438 (1995) 

[4] D. Boyanovsky, H.J. de Vega, R. Holman, D.S. Lee, and A. Singh, Phys. Rev. D 
51, 4419 (1995); D. Boyanovsky, M. D'Attanasio, H.J. de Vega, R. Holman, and 
D.-S. Lee, Phys. Rev. D 52, 6805 (1995). 

[5] S. Khlebnikov and I. Tkachev, Phys. Rev. Lett. 77, 219 (1996). 

[6] S. Khlebnikov and I. Tkachev, Phys. Lett. B390, 80 (1997). 

[7] S. Khlebnikov and I. Tkachev, [hep-ph/ 96 104771 , Phys. Rev. Lett. 79, 1607 (1997). 

[8] S. Khlebnikov and I. Tkachev, Phys. Rev. D 56, 653 (1997). 

[9] T. Prokopec and T. G. Roos, Phys. Rev. D 55, 3768 (1997). 

[10] L. Kofman, A. Linde, and A. A. Starobinsky, Phys. Rev. Lett. 76, 1011 (1996); 
I. I. Tkachev, Phys. Lett. B376, 35 (1996). 

[11] E. W. Kolb, A. D. Linde, and A. Riotto, Phys. Rev. Lett. 77, 3716 (1996). 
M. Yoshimura, |hep-ph/9605246| . 

[12] L. Kofman, A. Linde, and A. A. Starobinsky, [hep-ph/ 9704454 

[13] D. Polarski and A. A. Starobinsky, Class. Quant. Grav. 13, 377 (1996). 

[14] V. E. Zakharov, S. L. Musher, and A. M. Rubenchik, Phys. Rep. 129, 285 (1985). 



19 



